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Abstract 

The simple system composed of three neural-like noisy elements is considered. Two of them 
(sensory neurons or sensors) are stimulated by noise and periodic signals with different ratio of 
frequencies, and the third one (interneuron) receives the output of these two sensors and noise. 
We propose the analytical approach to analysis of Interspike Intervals (ISI) statistics of the spike 
train generated by the interneuron. The ISI distributions of the sensory neurons are considered to 
be known. The frequencies of the input sinusoidal signals are in ratios, which are usual for music. 
We show that in the case of small integer ratios (musical consonance) the input pair of sinusoids 
results in the ISI distribution appropriate for more regular output spike train than in a case of 
large integer ratios (musical dissonance) of input frequencies. These effects are explained from the 
viewpoint of the proposed theory. 

PACS numbers: 87.10.Ca, 87.19.1c, 43.72.Qr 



1 



I. INTRODUCTION 



Since 1980-th it is well known that noise in physical systems doesn't always play a negative 
role. The phenomena of coherence and stochastic resonance are found in different branches 
of science g and the typical Held of research regarding the constrnct,™ role of no.se is a 
wide class of neural systems. Indeed, it is very difficult to forget about noise, investigating 
various parts of the central or the peripheral nervous system, even if a study is carried 
out in a framework of some mathematical model: noise is inherent for the dynamics of the 
membrane potential (due to ion channels' noise) of neural fibers and somas; the synaptic 
junctions exhibit stochastic behavior; each neuron receives on the average 10 4 inputs from 
its neighbors [2j that in itself requires statistical methods of investigation, and etc. The 
motivating question of the given research is: how a signal survives in such noisy environment? 

Looking for an approach to the problem we concentrated our attention to sensory systems 
. Typically, in sensory systems there is a set of neurons (referred to as sensory neu- 
rons or sensors) receiving signals directly from the environment. For example, in a simple 
approximation of the mammals' auditory system it is supposed that the basilar membrane 
performs the Fourier transformation on an input sound signal 5|, and the sensory neurons 
attached directly to this membrane, percept different sinusoidal components (depending on 
coordinates of connection along the membrane) of the sound as input. Under driving of 
these signals and the mentioned noise they generate trains of short impulses (spikes), which 
are transmitted to other neurons (interneurons) along neural fibers. 

In a number of studies regarded to the neurodynamics under the noise influence the 
interspike intervals (ISI) statistics is of interest. In our model, composed of two sensors 
(stimulated by sinusoidal signals and noise) and one interneuron, we consider the ISI distri- 
bution (ISID) of each sensor to be known (from the previous works [6[]) and investigate ISIDs 
of the output signal of the interneuron driven by a mixture of noise and the sensors'spike 
trains weighted by coupling coefficients. The system with the similar structure has been in- 
vestigated in Ref. [§} for the purpose of Ghost Stochastic Resonance (GSR) phenomenon 
detection. The GSR term denotes existence of the maximum in the system response at some 
frequency, which is absent in the spectrum of the input signal. The maximum takes place 
at some optimal intensity of noise, which affects a system as well {9]. 

Though we have not investigated this phenomenon in the presented work, our topic is 
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closely connected with GSR studies due to high complexity (multimodality) of interneuron's 
ISIDs comprising peaks inappropriate to input sinusoids' periods or their multiples. 

We show how the input signal composed of two sinusoids is transformed by the proposed 
noisy system into different types of spike trains, depending on the ratio of input frequencies. 
Looking for the differences in the statistical sense, we find out that the output ISIDs for some 
combinations of frequencies have sharp shapes similar to ISIDs of an interneuron driven by 
a well recognized (on a noise background) regular signal. Also, there is another type of the 
output ISID (for the other frequencies combinations), which has a blurred shape similar to 
an ISID of a neuron driven by randomly distributed impulses. 

In fact, the difference between "sharp" and "blurred" shapes of ISIDs is more quantitative, 
than qualitative, but this difference indicates higher stability to noise of one combination of 
input sinusoids in comparison with another one. Investigation of this phenomenon can help 
to understand which types of input signals are able to survive in the noisy environment of 
the brain, which principles control this process, and what it means from the perceptional, 
cognitive, and other points of view. 

On the other hand, in the real life a human deals with relatively simple combinations of 
sinusoidal signals, when listens to music. It is well known that musical accords (combina- 
tions of tones) are classified as consonant (pleasant, harmonious) or dissonant (unpleasant, 
disharmonious), depending on the ratio between frequencies Thus, use of musical no- 
tations appears to be convenient in the context of our work for input signals classification 
purposes. However, we should emphasize that our results are obtained using the so-called 
"just intonation" musical accords, i.e. frequencies of input sinusoids are related by ratios 
of whole numbers, that is not appropriate for modern music, but is more suitable in the 
presentation sense. 

It is important also, that the consonance and the dissonance of accords are recognized 
by animals (which never deal with music) as well 4]. So, the underlying principles seem to 
be common and fundamental for the auditory neural system of mammals. This is the good 
reason to use the neural-like model of the auditory apparatus as the object of research into 
effects related with simple signals (like simple musical accords) propagation through a noisy 
nonlinear environment. 

It should be emphasized, that the "noise benefits" phenomena like coherence resonance, 
stochastic resonance, ghost stochastic resonance, etc. are appropriate candidates for a so- 
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lution of signal propagation and signal "survival" problems. But they allow to reveal a 
very particular peculiarities of signal propagation through the nonlinear noisy environment 
of neural-like systems and don't provide a full statistical picture. So, the main goal of the 
paper is to present an analytical description of principles, which control the statistics trans- 
formation process for spike trains propagated from one level of neurons to another one under 
the influence of noise. 

In the paper we first describe the chosen model in details. After that we propose the 
analytical description applied to the interneuron's ISI statistics. In order to prove the 
theoretical conclusions, we compare them with the results of computer simulations. Finally, 
we discuss an agreement of obtained results with the hypotheses of the consonance and the 
dissonance in music proposed by Helmholtz (1877) and Boomsliter&Creel (1961). 

II. MODEL 

As a basis for the investigated neural-like system we have chosen the widely used model 
called Leaky Integrate-and-Fire (LIF) neuron. The input neurons (sensors) are driven by 
the external sinusoidal signals, and the output one (interneuron) receives the weighted spikes 
of the input neurons. 

For simplicity we restrict consideration by a case of two sensors. As a result, the set of 
equations for our system can be written down in the following form: 



Here: v i(t) is the membrane potential of the i th sensory neuron; //j is the relaxation param- 
eter; Ai and f2; are the amplitude and the frequency of the corresponding harmonic input, 
respectively; Di is the sensor's noise intensity; (i = 1,2) are the independent sources 
of the zero-mean 5-correlated (< £i(t)£j(t') >= S(t — t')5ij) white Gaussian noise (WGN) of 
the sensors; v(t), fi, D, and £(£) are the membrane potential, relaxation parameter, noise 
intensity, and WGN of the output neuron (the third equation (jT|), < >= S(t — £'), 



The LIF neuron doesn't comprise any mechanism of spike generation. So, as soon as 
the membrane potential of any neuron of the model reaches the threshold value v t h we say 





(1) 



v 



< f (£) >= 0, and < £(t)&(tf) >= 0. 
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that the spike is generated at the threshold crossing instant of time. The corresponding 
membrane potential is reset simultaneously to the initial value: v for the sensors and v for 
the interneuron. 

The spike trains generated by the sensors and received by the interneuron are denoted as 

Ni(t) 

s i(t) = S ~ tij))i — 1)2. Each spike train is weighted by the corresponding coupling 

3=0 

coefficient Spikes are modelled by Dirac 5-functions. The instants of time ty correspond 
to threshold crossings by the sensors' membrane potentials, iVj(t) is the number of spikes 
generated by the i th sensor since the initial time. Obviously, the values t$j and their numbers 
Ni(t) are directly related with amplitudes of input signals that means the system is nonlinear 
(by the definition of nonlinearity) even though it is not clear from the model Eq.flTJ). 

All simulation and theoretical results presented in the paper are obtained using the fol- 
lowing set of constant parameters: fii — ^ = 1, A* — 0.3665, D\ — D% — D — 1.6 ■ 1CT 3 , 
ki = k 2 = 0.98, v i = v ° = 0, v° = —1, and v t h = 1, unless stated otherwise. 

For the output neuron the refractory period {T re f) is introduced explicitly: this neuron 
does not respond on any external signal after reset until the varying potential y°e~^ t ~ treaS > 
reaches the level v = —0.1. Hence, the refractory period can be written down in the following 
form: 

T ref = -\n(-10v°) . 
For the chosen parameters we have T re f = 6.28. 



III. THEORETICAL STUDY 



The first two equations of the system (JTJ) are, obviously, independent differential equa- 
tions, modelling the well-known Ornstein-Uhlenbeck process with the harmonic inhomo- 

geneity and the reset rule. The statistics of interspike intervals in this case can be obtained 

fl 

analytically or numerically [6j and we consider it to be known. 

The very important thing is that the spike trains of our sensors are non-Poisson ones. 
These spike trains are the input into the third neuron, and it means that the dynamics of 
the output neuron membrane potential is non-Markovian llfl . Hence, we are compelled to 
investigate the ISI statistics of the output neuron using another analytical approach. 
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A. Solutions and presuppositions 



It is possible to obtain analytical solutions for Vi(t) and v(t) 121 ] : 

Vi(t) 



u i(*0i) 7=t=7 COs{Qit 0i + 

V 



X 



v(t) = v(t )e-^-^ + Y: kiSitf) + y/D((t). 

Ni(t) 



2 

E 

i=i 



Here: Si(£) = ^ e ^ ti: >' is a sum of decaying impulses evoked by spikes of the i th sen- 

t 

sory neuron; Q(t) = f e-^-^i^dH is the colored Gaussian noise (Ornstein-Uhlenbeck 

to 

process) with the variance 

= (Cf{t)) = 2^" (1 — e~ 2flit ) and the probability distribution 



1 / s 2 



For the output neuron we have the same forms of c 2 (t), and u>f(s). to? and to are 
the reset (spike generation) instants of time for the sensors or the interneuron, respectively; 
(pi = arctan (jj^J ■ 

The temporal realizations of membrane potentials of neurons allow us to understand the 
conditions of spike generation by the output neuron and to establish connections between 
these events and input signals. 

In order to perform the following analysis we utilize three main presuppositions: 

1. The input harmonic signals are subthreshold to the sensors, i.e. the amplitude A, and 
the frequency fli are in such a relation, that the signal cos(fijt) is not able to evoke 
a spike of the i th sensor without noise (Di = 0). From the solutions ([2]) we obtain 

' ' < v th . (3) 



2. Only one spike can be generated at each period of the harmonic driving force. But, 
at the same time, the spiking on each period is the most probable situation, and it 
means the relatively (to high relaxation parameter /ij. Formally, the condition can 
be written down as: 

1 2tt 
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FIG. 1: The typical membrane potential realization t>j(i) and the ISI distribution of the sensory 
neuron (Oj = 0.6, A{ = 1.165). The highest probability of a spike after t = is near one period of 
external force (t = 10). The probability of firing after two, three, etc. periods decreases rapidly. 

3. Each of coupling coefficients hi is less than the threshold membrane potential value 
v t h- It means that any separate incoming spike evokes also a subthreshold impulse 
of the output neuron's membrane potential v(t), i.e. spike generation is impossible 
without noise. At the same time, the sum of two coefficients is greater than vth- 

h,2 < v th , 
k 1 + k 2 > v th . 

B. Probability distribution for the output neuron spike 

Let's make some theoretical estimations. Initially all three neurons of the system Eq. ([1]) 
are reset, i.e. t>i(0) = v®, v 2 (0) = v 2 , and v(0) = v°. Since the starting time is t — 0, we 
measure the first interspike period of the output neuron as the first passage time. The first 
passage time probability distributions (FPTPD) are considered to be known for the input 
neurons: pi(t) and p2{t), respectively. It means once time is started, we know necessary 
characteristics of a spikes sequence, coming from the 1 st and the 2 nd neurons to the 3 rd 
(output) one. Spikes of the sensory neurons have a highest probability to appear first time 
at maxima of harmonic driving force (Aj/ \/Vtf + 1) cos^tj + fa). They have a narrow 
probability distribution near each of these maxima, and the probability of skipping one, 
two, etc. periods decays exponentially (see the Fig. CD). 

So, we may analyze a probability dP^{t) to find the 3 rd neuron spike inside the short time 
interval [t, t + dt\. 
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For the chosen parameters there are only 4 probable situations in which the output (the 
3 rd one) neuron generates a spike: 

1. upon receiving a separate spike of the 1 st neuron; 

2. upon receiving a separate spike of the 2 nd neuron; 

3. upon receiving a I s * neuron's spike on a background of a 2 nd, s one; 

4. upon receiving a 2 nd neuron's spike on a background of a l st 's one. 

The "separate" spike means that at the time of its incoming the 3 rd neuron's membrane 
potential v(t) is driven only by noise, i.e. any previous perturbation over the noise level is 
relaxed. 

The "background" of some incoming spike means that this spike was not able to make 
fire the interneuron, but perturbed its membrane potential. This background decays expo- 
nentially with the decrement \i until it becomes hidden by noise. 

Four described situations exclude each other, so, we may take them as independent 
probabilistic hypotheses. The probability of each hypothesis realization is directly connected 
with probabilities of the 1 st and 2 nd neurons' spikes generation on the same time interval 
[t,t + dt], and before it. 

The probability of the 3 rd neuron spike generation in each of these 4 cases depends on 
the coupling coefficients /c 12 , the noise intensity D and the membrane threshold v t h- 

Hence, we obtain the first term of the contribution into the probability of the 3 rd neuron 
spike generation inside the time interval [t, t + dt] in the following form: 



where dP\(t) is the probability of the I s * neuron spike generation on the interval [t,t + dt]. 
Probl^ki + y/~D£(t) > v t h^ is the probability of the output neuron spike generation under 
the influence of the l si neuron's spike (the 1 st hypothesis probability). For example, if k\ is 
too small, then the 1 st neuron spike is not able in practice to make fire the output neuron. 

The same is applicable to the 2 nd neuron spikes influence, what provides us with the 
second term. 

For the third hypothesis let's imagine the 2 nd neuron spike comes to the 3 rd neuron and 
doesn't make it fire. In this case v(t) performs a short "jump" (its height is equal to k 2 ) and 
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decays exponentially towards zero. According to our presuppositions, during this decay we 
can expect only the 1 st neuron spike. And it has a real chance to make fire the 3 rd neuron. 
This "real chance" is equal to Prob j/^ + A^e - ^* - *') + VD((t) > v th }, where [t',f + dt] is 
a short interval of born of the previous incoming spike. It is obvious, that if the instant 
of time t' is too far from the current one t, then there is no effect given by the previous 
incoming spike, and in this case a spike at the current time t is named the "separate" one. 
So, we don't need to take into account all previous time t'. If the previous 2 nd neuron's spike 
doesn't evoke the 3 rd neuron spike, then the first one is totally forgotten by the interneuron, 
when 

i.e. when the noise amplitude becomes equal to the decayed impulse (not spike) evoked by 
the 2 nd neuron spike. By this way we obtain the meaningful period of time to integrate over: 

v \VdJ 

We also understand, that the whole situation is as seldom as high is the probability of the 
3 rd neuron firing under influence of a separate spike from the 2 nd neuron. It can be reflected 
by the factor: 

(l-Prob{k 2 + VD((t)>v th ^ . 
This way we obtain the next (third) term of dP^{t): 

dP 1 (t) (l - Prob [k 2 + VD((t) > v th yj x 

x / dP 2 {t')Prob \h + k 2 e-M-V + VD({t) > v th \ . 

t-T 2 L ) 

The opposite order of 2-spikes sequence (the A th hypothesis) contributes the term of the 
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same form with exchanged indexes 2 «-> 1. And the whole expression is: 

p 3 (t) = Pl (t)Prob {h + VD((t) > v th } + 



+p 2 {t)Prob [k 2 + VD((t) > v th } + 

+ Pl {t) (l - Prob {fe + v^DC(0 > «tft}) 

t 

x / p 2 (t')Prob{k 1 + k 2 e-i' ( - t - t "> + 
±VDC(t) > vth}dlf+ 



x 



t-T 2 



+p 2 (t) (l - Pro6 {fci + VDC(t) > v th j 

t 

x / p 1 (tf)Prob{kie-i* t -«) + k2+ 

t-Ti 

+VD((t)>v th \dt'. 



X 



(4) 



Here we switch our attention to the probability densities: pi(t) = dPi(t)/dt. And one must 
remember, that everything is valid only for t > T re f. 



C. Hypotheses' probabilities 

In order to make the expression (j4]) more clear, we should focus on the coefficients denoted 
as Prob{. . .}. The common representation of this factor is: 

Prob{v(t) >v t h}, 

i.e. the probability of the threshold crossing by the output neuron membrane potential. 

After expiration of the refractory period and before any incoming spike the output neuron 
membrane potential is equal to the Ornstein-Uhlenbeck process realization: v(t) = \^D((t). 

Once an external spike is received (e.g. the I s * neuron spike), v(t) performs an immediate 
jump to the value k\ + y/D(^(t). Obviously, due to the infinity of the derivation of this 
jump the probability of the 3 rd neuron spike depends only on a current value of the noise 
realization. That's why we may simply write the following: 

ProbU{t)>{vth-h)/VD\= J _wf(s)ds = 

(v th -ki)/VD 

= ±eric{y/E(v t h ~ h)} ■ 



10 



Here 



exp (— /zs 2 ) 



is the stationary probability distribution of the noise amplitude, erfc(x) is the complemen- 
tary error function. The stationary form is chosen, because the refractory period is long 
enough, and any "jump" of v(t) without spike generation does not reset the noise compo- 
nent. 

Using the same line of reasoning, it is easy to understand, that in a case of 2 incoming 
spikes close in time we obtain almost the same result (e.g. the 1 st neuron spike comes on 
the background of the decaying "jump" evoked by the 2 nd neuron spike): 



When the previous incoming spike and the current one are close in time to each other, 
the difference (t — t') is very small. It can be almost equal to zero (simultaneous spikes). In 



On the other hand, when (t — t') is large (long period between the previous and the 
current incoming spikes), we find, that: 



and this is the minimum of considering probability as a function of difference t — t'. 





such a case we deal with the maximum of Prob < ((t) > (yth — k\ — A^e ^ *'-')/ V^D 




Denoting 




(5) 



and 
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FIG. 2: Interaction functions introduced in Eq. ([5]) for weak and strong connections: a) - weak 
connections: k\ = 0.6, k 2 = 0.7. There is a rather sharp boundary of a time period, when the 
second incoming spike is able to finish the "work" of the previous one. The separate sensor's spike 
is not practically able to fire the output neuron; b) - strong connections: k\ = 0.94, k 2 = 0.98. It 
is easy to see that any separate spike is able to make fire the output neuron. 

we can rewrite Eq. (jlj) as: 

p 3 (t) = pi(t)<f>oi(h) + P2(t)% 2 (k 2 ) + 



+pi{t) (1 - $02(^2)) / p 2 (t')$ 1 (k u k 2 ,t-t')dt'+ 

t-T 2 



(6) 



where 



+p 2 {t) (1 - $oi(*i)) / p l {tf)^ 2 (k 2> ki,t-1f)dtf, 



rp 1, fk h2 

T 12 = - In -^= 



The $1,2(^1,2, k 2 i,t — t') factors are depicted at the Fig. [2] as functions of the time difference 
t-t'. 



D. One more important multiplier 

Regardless of shapes of pi, 2 {t) the FPTPD ps(t) must have one important characteristics: 
if the output neuron spike appears at some earlier time, then this circumstance decreases the 
probability of the spike in all later moments. We can reflect such a property by multiplying 

12 



ps{t) by 



t 



1 - / P3 (t')dt' 



o 



The problem is that our p 3 (t) may occur to be not normalized after the previous calculations 
(we explain it below, in the section HVBj) . But it is possible (without losing of generality 



and facing with any contradictions) first to obtain p 3 (t) as described above, then normalize 
it, and then multiply the result by the mentioned multiplier: 



p 3 (t) = p 3 (t) \1-J p 3 (t')dt' | , (7) 
where p 3 {t) - is the normalized probability distribution. 



E. A set of the interneuron states 



Now we should recall that all previous calculations are valid until the first spike is gen- 
erated by the 3 rd neuron. The question is: what happens after? 

At the moment of the 3 rd neuron's spike generation its membrane potential v(t) is reset 
to initial value v°, and the interneuron "forgets" all previous history. We suppose that 
the 3 rd neuron spike is evoked by a spike of the 1 st or the 2 nd sensor exactly at the same 
moment. Let it be the 1 st sensor, which makes fire the output neuron. It is also reset to its 
initial membrane potential value v®. Consequently, after reset FPTPD pi(t) has the same 
shape as it was previously. The other sensor is not reset simultaneously with t>i(£) and v(t). 
Therefore, its FPTPD p2{t) is shifted now in comparison with the initial situation. The Eq. 

is valid, but now it provides us with new FPTPD p { 3 \t), where the index (.)(°) is used 
for the initial situation. The time is measured now since the moment of last spike generation 
by the interneuron. The same is correct after each reset of the interneuron: one of pi,2(t) is 
similar to its initial form, while another one is shifted to the left or right. Hereafter we say 
that the output neuron gets into some state after each reset. These states are defined by 
corresponding shifts of FPTPDs pi^if) from the "viewpoint" of the Eq. ([6]). For detailed 
description see the section IIVBI In the case of sinusoidal inputs and a finite number of 
sensors we have a finite number of these states. Hence, the resulting FPTPD of the output 
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neuron should be written as 

Pout{t) = a pf\t) + aipj s 1) (t) + 

+a 2 pP(t) + --- + a M . 1 pi M - 1 \t), 
where M is a whole number of the interneuron's states. The coefficients denote relative 
frequencies of switching into according states. We talk about values of these coefficients 
below in the section TlV B[ Since the distribution p^\t) is particular for the k th state, it 
would be incorrect to use here p^\t) (without a hat) and then to perform the same operation 
as in Eq. t^. 



F. Example for different frequencies of sinusoidal inputs 

Let's take two sensors with input sinusoidal signals of different frequencies fli 7^ We 
only suppose, that these frequencies are in a ratio of some integers m and n, i.e. Q1/Q2 — 
m/n. This means that the first (m — 1) peaks of pi(t) don't coincide with the first (n — 1) 
peaks of P2(t). And the m th peak of pi{t) coincides with the n th peak of P2{t). 

Consequently, the output neuron has M = (m — 1) + (n — 1) + 1 = m + n — 1 different 
possible states against peaks of pi^t), i.e. if it is reset together with any spike of the first or 
the second sensor, then since the reset time there is always only one of M different variants 
of an incoming spike train (the superposition of spike sequences from both sensors) with a 
definite probability density in time for each incoming spike. 



IV. NUMERICAL EXPERIMENTS 

A. Consonance and Dissonance in music 



The sinusoidal signal is considered to be the simplest one in an investigation of different 
systems. May be the first level of complication is the linear combination of two sinusoidal 
signals of different frequencies. And here we face with the set of very old questions related 
with musical accords. 

The Pythagoreans discovered that the accord of two sinusoidal signals sounds pleasant 
(consonant) if their frequencies ratio is m/n, where m and n are the small integers (e.g. 2/1, 
3/2, 4/3) 



10j. Conversely, if m and n are the large numbers (e.g. 45/32), then the accord 



sounds dissonant, i.e. unpleasant. 
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FIG. 3: The consonant accords: under each picture there is the ratio of frequencies (m/n) and the 
name of the accord used in the common musical terminology. All curves are obtained through the 
direct numerical simulation of the system Eq. (pQ) with 0,2 = 0.6, A2 = 1.165, f2i = (rn/n)^, and 
A\ according the subthreshold input sinusoidal signal condition Eq. ([3|). 



In the context of our investigation it is really interesting, how the dissonance and the 
consonance are mapped to ISI distributions. 

In the Figs. [3] and H] there are the distributions for consonant and dissonant accords, 
respectively. It is easy to see the higher integers m, n the regularity less in an appropriate 
distribution of ISI, although the structure of the input signal is in principal the same: two 
sinusoids. These curves are obtained through the direct numerical simulations of the system 
Eq. (TjQ). The theoretical part of our work is focused on building a basis for the simulations 
results. 



B. Verification of theoretical conclusions 

The formula (j5D is obtained under a set of assumptions. So, this theoretical result should 
be compared with the results of numerical experiments. Here we present the idea of usage 
of the expression and check its validity. 

Let's take, for instance, the "Perfect 4 th " accord, which consists of 2 sinusoids of frequen- 
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FIG. 4: The dissonant accords: under each picture there is the ratio of frequencies (jn/n) and the 
name of the accord used in the common musical terminology. All curves are obtained through the 
direct numerical simulation of the system Eq. (fTJ) with 0,2 = 0.6, A<i = 1.165, f^i = (m/n)02, and 
A\ specified by the subthreshold input sinusoidal signal condition Eq. (J3j) - 



cies related by the ratio Q1/Q2 = 4/3. 

The FPTPDs pi(t) and p2{t) are known to us from the numerical simulations of sensors, 
for example. If the figures of these distributions are placed in a column (Fig. [5J State 0), 
then it is easy to see, that the 3 rd peak of P2(t) coincides with the 4 th peak of pi(t). All other 
peaks don't coincide, and (as it is explained in the section IHI Fj) here we have 4 + 3 — 1 = 6 
different possible states of the 3 rd neuron: State 0, State 1, . . . , and State 5. Let's establish 
the correspondence between these states and the peaks of pi,2(t) as it is shown on the Fig. 
[5j numbers of the states are placed into circles. Area under each peak means the probability 
to find an incoming spike at the defined short period of time. If this spike evokes the spike 
of the 3 rd neuron, it is switched into the appropriate state. 

Initially all three neurons are reset. This is the th state. The most probable and close 
in time spike comes from the 1 st sensor (Fig. [5]). If this spike makes fire the interneuron, 
then it is switched into the I s * state, where the most probable and close spike comes from 
the 2 nd sensor. Obviously, this spike (if it is generated) comes during the refractory period, 
so, the closest valuable spike in the I s * state comes again from the 1 st sensor and has the 
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FIG. 5: FPTPDs of three neurons in the states and 1. This is the example of the probable 
transition between two states for the case of input sinusoids with frequencies related by ratio 4/3 
(the Perfect 4 th accord). P2(t) is shifted in the State 1 in comparison with P2(t) in the State 0. As 
a result, pg (t) and p% (t) are different. All possible states of the interneuron are: State 0, State 
1, . . . , and State 5. The peaks of pi,2(i) are marked by numbers in circles in order to establish the 
correspondence between them and the interneuron states. 

possibility to switch the interneuron into the 3 rd state, and etc. 

Here we recall the section [III Dl and notice that the peak of pi(t) or P2(t), which appears 
during the refractory period, is the main reason why we may find p%(t) to be not normalized 
in Eq. (j3J). This "invisible" to the output neuron peak does not contribute into the ps(t) 
peaks, but it is the big meaningful part of an appropriate normalized sensor's FPTPD. 

The analysis of p$(t) peaks in each state shows that in the case of musical accords and 
strong connections the interneuron gets into all possible m + n — 1 states almost uniformly, 
i.e. all states make almost equal contributions into the common p ou t(t). 
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FIG. 6: ISI distributions of the output neuron for different accords. Solid line is the theoretical 
result. Dashed line is the distribution obtained throughout the direct numerical simulation of the 
system ([TJ. The parameters are /ii = fj, 2 = 1> A 4 = 0.3665, k\ = k 2 = 0.97, D\ = D 2 = D = 
1.6 • 10~ 3 . Perfect 4th: A 1 = 1.165, fti = 0.6, A 2 = 1.085, fi 2 = 0.45. Minor 3rd: A l = 1.125, 
fti = 0.54, A 2 = 1.085, n 2 = 0.45. Major 2nd: A 1 = 1.2, Qt = 0.675, A 2 = 1.165, VL 2 = 0.6. 

Therefore, the simplest way to obtain the final output distribution is to directly sum all 
Pz\t) and then to normalize this result. In other words, all coefficients (see the section 
IIII Ej) can be set to unit: 

M 

Pontic = ^ ■ (8) 

k=l 

These approximate conclusions provide p ou t{t) curves very similar to ones obtained 
through direct numerical calculations of the system ([ID. The examples of compared results 
are shown in the Fig. [6j 



18 



V. THE CALCULATION ALGORITHM 



Summarizing the previous sections, let's present the described theoretical approach in 
the form of the calculation algorithm. Thus, in order to obtain the interneuron's ISID curve 
under chosen parameters of the system Eq. ([T]) we should perform the following steps. 

1. To obtain the sensors' FPTPDs pipit) using the direct numerical simulation of the 

r *.— — . — . — — 

2. To find all possible States of the interneuron against the peaks of pi,z(t) (see Fig. [5] 
for example). 

3. To calculate p^\t) in each State using Eqs. (0[7]). 

4. And, finally, to sum and normalize the calculated pij fit) in accordance with Eq. ([5]). 

Despite of the relative complexity of the algorithm, its usage decreases consumption of 
resources necessary for smooth interneuron's ISID obtaining. It allows also to perform fast 
estimations and provides anyhow the consistent theoretical description of the noisy nonlinear 
system Eq. ([!]). 



VI. HYPOTHESES OF CONSONANCE AND DISSONANCE 

There is a few of main hypotheses explaining why animals, including humans, feel har- 
mony or disharmony listening to different tones combinations. We suppose that the input 
signals, which are transformed into spike trains with blurred, i.e. noise-like Interspike In- 
tervals Distributions, are felt unpleasant (dissonant, inharmonious) due to the analysis, 
recognition and survival in noisy environment of the brain problems. 

Let's try to understand the correlation between that viewpoint and some other hypotheses 
of the dissonance. 

n 

Helmholtz (1877) [5| proposed the notion that dissonance arises due to beating between 
adjacent harmonics of complex tones. In effect, dissonance arises due to rapid amplitude 
fluctuations. 

It is possible to prove mathematically that, if the input frequencies into our system are 
in ratio Q1/Q2 = m/n (where m > n), then the minimal distance between peaks of p^\t) 
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is T min = (27r/f2 2 )/m = T 2 /m (see the Fig. [5j the peaks 1, 2 of the p% (t), and the peaks 
4, 5 of the Ps(t)) that defines the distance between peaks of the final ISID p ou t(t)- That's 
why the sufficiently high value of m means the blurred ISID, typical for dissonant accords. 
In such a way, we show that even for pure input tones, and even if they are not close in 
frequency in order to produce beats, we may feel the dissonance. Hence, the hypothesis by 
Helmholtz continues to be correct, if we look at the minimal distance among all peaks of 
Pi(t) and P2{t), and not only the firsts ones, which show the difference between the input 
tones' frequencies. 

Another theory is the Long Pattern Hypothesis from Boomsliter and Creel (1961) [3] 
which states that a consonance is based on the length of the overall period of a stimu- 
lus. They show that consonant intervals, based on simple integer ratios of fundamental 
frequencies, have shorter overall periods than do dissonant intervals. 

Indeed, as we obtain for our model, the higher integers are m and n, the higher number 
of states (m + n — 1) the interneuron has against the pattern of pipit) peaks. In fact, the 
sequence of the states repeats periodically in time with the period T sta te = (2ir/Qi)m = 
(2-7r/f22)n, which is the period of phases coincidences of cos fl\t and cosf^t, i.e. the overall 
period. But the interneuron gets into each state randomly. So, for the high number of 
states (dissonance) it is necessary much time in order to recognize some regularity inherent 
to the output spike train. Conversely, in the case of consonant input, the same amount of 
the spike train statistics details can be acquired in shorter periods of observation. Thus, the 
consonant inputs are in the priority against the dissonant ones, from the analysis viewpoint. 



VII. DISCUSSION AND CONCLUSIONS 



In the work we try to follow a signal propagating throughout the neural-like system. The 
second layer of the system doesn't allow applying the framework of Markovian processes. 
Nevertheless, we propose the qualitative analysis yielding the main result of the work: the 
analytical expressions and the consistent algorithm applicable for an investigation of the ISI 
statistics and its transformations. 

The proposed algorithm is ready to be used for quick estimations of output distributions 
because of step-like shapes of the functions called "$(. . . )" and narrow peaks of FPTPDs 
pi(t) and p 2 (t). 
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On the other hand, the found procedure is clear enough to be implemented in the widely 
used programming environments. In such an implementation it provides a rather precise 
approximation (see the Fig. E]) of the output ISI distribution (ISID) given by the direct 
computer simulation of Eq. ([T]). We should also emphasize that in order to obtain a smooth 
curve of the ISID, using the direct numerical simulation of the system ([Tj), it is necessary to 
consume much more temporal, soft- and hard- ware resources than in the case of a program 
implementation of the proposed algorithm usage. 

In the simple case of the auditory system model we are able to discover existence of some 
accords (a combination of two sinusoidal signals), which evoke ISIDs blurring very fast with 
propagation from one neural layer to another (Fig. H]). And in our study these accords are 
the same as the dissonant ones in music, i.e. the dissonant accords are the ones, which are 
not able to "survive" in the noisy neural environment after a number of interneurons layers. 

We also show that from the perceptional point of view the dissonant accord's ISI statis- 
tics needs more time to be collected in comparison with the consonant accord's one. The 
latter one evokes a sharp ISID's shape, which is able to "survive" a number of proposed 
transformations, i.e. the same algorithm is applicable in order to understand what happens 
to the consonant accords on deeper layers of the neural system. 

As it is easy to see the output ISID contains peaks corresponding to quasi-periodical 
spike generation at frequencies, which are absent in the input signal. So, it is possible and 
intersting to investigate the Ghost Stochastic Resonance phenomenon [7, 9] in details for this 
model. However, the current paper is focused on the theoretical approach to the whole ISID 
picture shaping. All sophisticated tuning of coupling coefficients, input frequencies, and 
noise intensities can be performed separately in a sake of resonances investigation, and this 
analysis can be also augmented by results revealed from the proposed analytical approach. 

The obtained results may be applied also in the context of such recent studies as, for 
example, the stimulus reconstruction from neural spike trains, where the information trans- 
mission under the noise influence is investigated [14]. The other suitable context of these 
results application is the continuous investigation of the neuron's behavior under the in- 
fluence of a constant bombardment of inhibitory and excitatory postsynaptic potentials 
somehow resembling a background noise that is typical for functioning conditions of, for 
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example, the neocortical neurons 
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